Global firing induced by noise or diversity in excitable media 
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We develop a theory for the emergence of global firings in non-identical excitable systems subject 
to noise. Three different dynamical regimes arise: sub-threshold motion, where all elements remain 
confined near the fixed point; coherent pulsations, where a macroscopic fraction fire simultaneously; 
and incoherent pulsations, where units fire in a disordered fashion. We also show that the mechanism 
for global firing is generic: it arises from degradation of entrainment originated either by noise or 
by diversity. 

PACS numbers: 05.45.-a, 05.45.Xt, 05.40.-a, 02.50.-r 
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Excitable behavior appears in a large variety of phys- 
ical, chemical and biological systems [l|, |2j- Typically 
this behavior occurs for parameter values close to an os- 
cillation bifurcation, and is characterized by a nonlinear 
response to perturbations of a stationary state: while 
small perturbations induce a smooth return to the fixed 
point, perturbations exceeding a given threshold induce 
a return through a large phase space excursion (firing), 
largely independent of the magnitude of the perturba- 
tion. Furthermore, after one firing the system cannot be 
excited again within a refractory period of time. In many 
situations of interest, the firings are induced by random 
perturbations or noise 

In coupled excitable systems, global firing (a macro- 
scopic fraction of the units fires simultaneously) excited 
by noise has been observed in chemical excitable me- 
dia H, neuron dynamics and electronic systems @, 
and it has been described through several theoretical ap- 
proaches 0, ■ This synchronized firing can be consid- 
ered as a constructive effect induced by the noise. Other 
examples in which noise actually helps to obtain a more 
ordered behavior are stochastic resonance 0, stochas- 
tic coherence (or coherence resonance) lljj, an d noise- 
induced phase transitions flll |. 

Diversity, the fact that not all units are identical, is 
an important ingredient in realistic modeling of coupled 
systems. Ensembles of coupled oscillators with dive rsity 
have been paradigmized and largely studied 0, Il4| . 
with the result that synchronized behavior can appear 
once the disorder induced by the diversity is overcome by 
the entraining effect of the coupling. It has been shown 
that in a purely deterministic excitable system diversity 
may induce global collective firing 0] if a fraction of 
the elements are above the oscillatory bifurcation. So, 
diversity and noise might be expected to play a similar 
role. 

In this work, we develop an analytical understanding 
for the emergence of global firing in coupled excitable 
systems in presence of disorder, either noise or diversity. 
We show that three different dynamical regimes are pos- 
sible: sub-threshold motion, where all elements remain 



confined near the fixed point; coherent pulsations, where 
a macroscopic fraction fire simultaneously; and incoher- 
ent pulsations, where units fire in a disordered fashion. 
Remarkably, the coherent behavior appears through a 
genuine phase transition when the noise intensity, the 
coupling or the diversity cross a critical value. A second 
phase transition to the disordered (incoherent) phase is 
recovered for large enough noise intensity or diversity, or 
small enough coupling. The mechanism for global firing 
is the degradation of entrainment which can be originated 
either by noise or diversity. This is generic and opens a 
new scenario for experimental observations. 

We consider as a prototypical model an ensemble of 
globally coupled active-rotators <f>j{t), j — 1, ...,JV, 
whose dynamics is given by [lri | 
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The natural frequencies ujj are distributed according to 
a probability density function g(u>j), with mean value u 
and variance a 2 . Notice that u>j < 1 (resp. u>j > 1) corre- 
sponds to an excitable (resp. oscillatory ^t|) behavior of 
the solitary rotator j. Throughout the paper we consider 
the case u> < 1. D is the intensity of the Gaussian noises 
£j of zero mean and correlations (£j = 8{t—t')8jk, 

and C is the coupling intensity. 

To characterize collective behavior we use the time- 
dependent global amplitude, p{t), and phase, H>(t) 0, 
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The Kuramoto order parameter p = (p(t)), where (•) de- 
notes the time average, is known to be a good measure 
of collective synchronization in coupled oscillators sys- 
tems, i.e. p = 1 when oscillators synchronize <j>j(t) = 
<fik(t), Vj, k, and p — > for desynchronized behavior. No- 
tice, however, that the Kuramoto parameter adopts a 
non-zero value even when all the variables <f>j , being equal 
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to each other, are at rest. To discriminate between this 
static entrainment from the dynamic entrainment (time 
synchronization) of excitable systems when all units fire 
synchronously, we use the order parameter introduced by 
Shinomoto and Kuramoto 0] 
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which differs from zero only in the case of synchronous 
firing. Finally, a measure for the activity of the units, 
widely used in problems of stochastic transport in non- 
symmetric potentials is the current 
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A non-zero current J describes a situation in which the 
systems are firing (not necessarily synchronized). 

We now provide an analytical theory to understand the 
behavior of p, £ and J as a function of the control param- 
eters, C, D and a. The theory proceeds in three steps. 
First, under the assumption of entrainment, we derive a 
dynamical equation for the global phase 'J, depending 
on the value of the Kuramoto parameter p. Second, us- 
ing the solution of that equation, we obtain expressions 
for £ and J which depend on p. Finally, we calculate 
self-consistently the value of p. 

Averaging equation Q over the whole ensemble and 
using the definition of global amplitude and phase of 
Eq. ((2J) we have 
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where £(t) is a Gaussian noise of zero mean and correla- 
tions (£(£)£(£')) = 6(t — t'). Taking the time-derivative of 
Eq. J5J) and introducing Sj(t) = 4>j(t) — we obtain: 
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We consider now that the rotators are synchronized in the 
sense that 6j (t) <C 1 and substitute the expansion e l5k = 
1 + iSk + 0(<5fe) in the previous expression. Equating real 
and imaginary parts, we obtain 
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The definition of 6i leads to p(t) = N^ 1 £\ 
p(t) — 0(6%) and, consistently with the order of the ap- 
proximation, we can replace in equation Q the time 
dependent p(t) by the constant value p. Therefore, 
Eq. |JSJ| can be approximated by p^(t) = uj — psm^(t) + 
y/D/N£(t), which in the limit N — > oo, reduces to 



f(t) = sin*(t). 
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It is remarkable that the global phase obeys the same 
dynamics than the individual units but with a natural 
frequency scaled with p, the Kuramoto parameter mea- 
suring the entrainment degree. Therefore, a decrease in 
the entrainment lowers the global excitability threshold 
from uj — 1 to u> — p and the system can start firing syn- 
chronously. The effect can be understood as a broadening 
of the distribution of the phases cf>, so that a fraction of 
the rotators crosses over the threshold and, if the cou- 
pling is large enough, they pull a macroscopic fraction of 
the oscillators. Thus degradation of the entrainment has 
the paradoxical effect of increasing the coherent firing. It 
is essential to realize that Eq. (JSJ depends only on the 
value of p and not in the specific way the degradation of 
p is achieved, so that similar effects can be achieved ei- 
ther increasing the noise, either decreasing the coupling, 
or increasing the diversity in the natural frequencies; a 
significantly insightful result not previously understood 
nor discussed. 

Let us now express £ and J as a function of the Ku- 
ramoto order parameter p. In the case p < uj, the solution 
of © is given by 
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where f2 = yj (uj/ p) 2 — I is the frequency of the global 
phase oscillations. The current is obtained from Eq. JSJ, 
J = uj — (p sin(\&)). Time averages are computed over a 
period T = 2n/Q using Eq. (JS) 
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For p > uj, J — 0. 

Approximating p(t) by a constant value, 
the Shinomoto-Kuramoto parameter £ = 
p {le^w — (e 1 *^)]) can be computed performing 
again the time averages over a period T using Eq. tfify: 

c = - V 2 ( w - V^7){" + p)k (^-] , (ii) 

7T V \p — UJ J 

where K(m) is the complete elliptic integral of the first 
kind H3. If p > uj we get C = 0. 

As a final step, we derive a equation for p using a self- 
consistent, Weiss-like, mean field approximation, which 
assumes constant values for the global magnitudes and 
then averages over their probability distribution nam. 
For our particular case, we start by rewriting Eq. QJ as 
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where we have defined the potential 

V ((/); p, uj) = -uj(j) - cos(0) - Cp cos(* 



(13) 
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FIG. 1: Symbols represent p, £ and J as obtained numerically from Eqs. Q. Solid lines are the theoretical results. Panel 
(a) shows the variation with respect to the noise intensity D in absence of diversity, a = 0, for frequency lo = 0.95, coupling 
C = 4, and different system sizes: N = 50 (o), N = 10 2 (x), N = 10 3 (A), N = 10 4 (o). Panel (b) displays the same results 
as a function of C for D — 1.0. Panel (c) shows the variation with respect to diversity a for D = 0.3, C — 4 and being a 
uniform distribution. In all cases there are three regimes: (I) no firing, (II) synchronized firing and (III) unsynchronized firing. 



Note that the coupling appears only through the global 
parameters p and For fixed p and the stationary 
probability distribution function reads [21j 

Jo 
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where Z is a normalizing constant. From its definition, 
we have p = (l/N) J2k=i ( cos (0fc — an d we obtain 

P = 
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where we have performed a triple average: with respect 
to the distribution JTIJl, with respect to the distribution 
g(uj) of natural frequencies and with respect to the dis- 
tribution P(^>;p) of the global phase which is inversely 
proportional to the instantaneous velocity given by the 
dynamics (JSJ, namely P( X F; p) = (1 / '2ir) \J 'lu 2 — p 2 /(ui — 
psin*) for p < lu and P(^>;p) = — arcsin(Lj/p)) oth- 
erwise. The self-consistent equation (|15f) for p needs to 
be solved numerically. 

In the following, we discuss the theoretical results and 
compare them with the numerical results obtained from 
a numerical integration of Eqs. Q). Fig. QJi shows p, 
£ and J as function of the noise intensity D in absence 
of diversity. The solid lines correspond to the theoreti- 
cal results while symbols show the numerical results for 
different system sizes. In this figure, we can observe the 



three aforementioned behaviors: For small noise intensity 
(regime I) each rotator fluctuates around its fixed point. 
Although for un-coupled rotators noise would eventually 
excite some spontaneous random firings, the coupling of a 
large number of units suppresses these individual firings. 
The Kuramoto parameter p close to 1 and the deviations 
from unity are due to the small dispersion induced by 
noise. Region I is, in fact, characterized by p > lu for 
which our theory predicts that the Shinomoto-Kuramoto 
parameter £ and the current J vanish which physically 
reflects the nonexistence of collective movement. In this 
region, the numerical results for p, £ and J are in excel- 
lent agreement with the theoretical predictions. 

Our theory predicts that a transition to a dynami- 
cal state characterized by synchronized firing behavior 
(regime II) takes place when p — lu, in very good agree- 
ment with the numerical results. This transition is clearly 
signaled by non- vanishing values of £ and J. The predic- 
tion of p is good for a large part of region II (up to values 
of p = 0.7). Later it underestimates its value. 

For very large noise intensity, the rotators desynchro- 
nize while keeping a non-zero current value (regime III). 
Hence, the synchronized activity, as measured by £ goes 
though a maximum as noise amplitude increases. Our 
theory predicts that the transition between regions II 
and III occurs for ,0 = where £ = and the current 
takes the maximum possible value J = lu. Surprisingly, 
since the small values of p in this transition are beyond 
the assumptions of the theory, the location of the second 
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transition is also well predicted. Moreover, the whole 
shape of the Shinomoto-Kuramoto parameter £ is well 
reproduced over the whole range. The maximum of Eq. 

occurs for p sa 0.821a;, which is well confirmed by the 
numerical results. The theoretically predicted current J 
fits the numerical values in the same range than p. Note, 
however, the numerical simulations show a local maxi- 
mum for the current J which indicates a local increasing 
in the total transport due to the coherent dynamics in 
the regime II. This local maximum is not present in the 
theoretical approximation. 

Some of these states were already described by Ku- 
ramoto and Shinomoto ^^|- By looking at the probabil- 
ity distribution of <f>i, these authors identify two regions 
in parameter space: the time-periodic regime (P) and 
the stationary regime (S). Region P corresponds to our 
regime II where the order parameter ( is different from 
zero and there is collective motion of the oscillators. Our 
findings allow us to split region S of these authors in our 
distinct regions I and III: while region I is a fluctuat- 
ing regime around the steady state, region III has a high 
activity as characterized by a non-zero current J. 

These results indicate that noise acts in two antago- 
nistic ways: while a given noise intensity can excite the 
sub-threshold units, forcing a synchronized firing, large 
amplitude noise deteriorates the synchronization proper- 
ties of the ensemble. This scenario resembles the so called 
noise induced phase transitions [TlT ] in which a transition 
to an ordered ferromagnetic-like state is induced by in- 
creasing the noise intensity; the order is destroyed again 
for large enough noise. Here, the transition is towards an 
organized collective motion of the active rotators. 

The reverse scenario can be observed varying the cou- 
pling strength C, see Fig. The Kuramoto parame- 
ter p increases with C (notice, however, the existence of 
a small bump in the numerical results), indicating that 
the degree of synchronization increases with coupling, as 
expected. A large coupling suppresses noise-induced fir- 
ings, and the system is macroscopically at rest, regime 
I, as indicated by the vanishing of £ and J. For weak 
coupling the noise induces desynchronized individual fir- 
ings (regime III) characterized by a macroscopic current 
J and again a zero value for For intermediate values of 
the coupling (regime II) the interplay of noise and cou- 
pling leads to the largest degree of synchronized firing 
with a large value for 

Finally Fig^i shows p, C and J as a function of the 
diversity a. It is clear in the figure the existence of the 
same three regimes that were obtained by varying the 
noise intensity or the coupling. Altogether Fig. ^ dearly 
illustrates the fact that similar effects can be achieved 
increasing the noise, decreasing the coupling or increasing 
the diversity in the natural frequencies as theoretically 
predicted. 

In summary we have developed a theory for the emer- 
gence of coherent firing in subthreshold excitable sys- 



tems. It arises as a consequence of a phase transition 
induced by noise or diversity. The underground mech- 
anism is the degradation of entrainment originated by 
disorder. This leads to the somehow paradoxical result 
of establishing a lower effective threshold for collective 
firing. It does not matter the specific source of disorder 
responsible for entrainment degradation, either noise or 
diversity when nonlinearly combined with coupling lead 
to similar results. Coherent firing is thus achieved by a 
global saddle-node bifurcation, triggered by the degree 
of entrainment of the ensemble, p. This mechanism is 
not restricted to the model we considered, it will exist 
in any physical, chemical or biological system with the 
necessary generic ingredients. 
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